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ABSTRACT 

We present the first self-consistent model for the dynamics and the radiative transfer occurring in 
bright X-ray pulsar accretion columns, with a special focus on the role of the shock in energizing the 
emerging X-rays. The pressure inside the accretion column of a luminous X-ray pulsar is dominated 
by the photons, and consequently the equations describing the coupled radiative-dynamical structure 
must be solved simultaneously. Spectral formation in these sources is therefore a complex, nonlinear 
phenomenon. We obtain the analytical solution for the Green's function describing the upscattering 
of monochromatic radiation injected into the column from the thermal mound located near the base of 
the flow. The Green's function is convolved with a Planck distribution to model the X-ray spectrum 
resulting from the reprocessing of blackbody photons produced in the thermal mound. These photons 
diffuse through the infalling gas and eventually escape out the walls of the column, forming the 
observed X-ray spectrum. We show that the resulting column-integrated, phase-averaged spectrum 
has a power-law shape at high energies and a blackbody shape at low energies, in agreement with the 
observational data for many X-ray pulsars. 

Subject headings: pulsars: general — stars: neutron — shock waves — radiation mechanisms: non- 
thermal — methods: analytical — X-rays: stars 



1. INTRODUCTION 

The first pulsating X-ray sources, Her X-l and Cen 
X-3, were discovered in 1971. Since then, over 50 
sources have been detected in the Galaxy and the Mag- 
ellanic Clouds, with luminosities in the range Lx ~ 



10 



34-38 



ergs s and pulsation periods in the range 



0.1 s < P < 10 3 s. In many cases the X-ray spectra are 
well fit by a combination of a power-law spectrum plus a 
blackbody component (e.g., Coburn et al. 2002; di Salvo 
et al. 1998; White et al. 1983). Most spectra also dis- 
play quasi-exponential cutoffs at E ~ 20 — 30 keV, and 
in some cases cyclotron lines are observed near the cut- 
off energy, implying magnetic field strengths ~ 10 12 G. 
Iron emission lines also appear to be present in certain 
sources. The X-ray emission from bright X-ray pulsars 
is powered by accretion of material from the "normal" 
binary companion onto the neutron star. 

Efforts to calculate theoretical X-ray pulsar spectra 
based on cither static or dynamic models have gener- 
ally yielded results that do not agree very well with 
the observed power-law profiles (e.g., Meszaros & Nagel 
1985a,b; Nagel 1981; Yahel 1980; Klein et al. 1996), 
which is probably due to the neglect of the central role 
of the shock in upscattering the radiation. Consequently 
no satisfactory theoretical model for the production of X- 
ray pulsar spectra has yet emerged (Coburn et al. 2002). 

Motivated by the problems with the current models, 
we reconsider here the physical picture originally pro- 
posed by Davidson (1973), in which the accreting gas 

1 pbecker@gmu.edu 

2 also Department of Physics and Astronomy, George Mason 
University, Fairfax, VA 22030-4444, USA 

3 michael.wolff@nrl.navy.mil 



passes through a radiative, radiation-dominated shock 
before settling onto the surface of the star. In this situa- 
tion, most of the photons are produced in a dense "ther- 
mal mound" at the base of the column, just above the 
stellar surface. These relatively low-energy photons are 
upscattered in the shock and eventually diffuse through 
the walls of the column, carrying away the kinetic en- 
ergy of the gas. This Fermi mechanism characteristically 
produces a power-law continuum at high energies. 

In this Letter, we demonstrate for the first time how 
the accretion shock itself can directly produce spectra 
very similar to the continuum emission seen in many X- 
ray pulsars. Although this approach has been suggested 
before (e.g., Burnard, Arons, & Klein 1991), it has not 
been carried out quantitatively. The work presented here 
therefore represents a significant step in the development 
of a detailed theory for the spectral formation process in 
X-ray pulsars. 

2. DYNAMICS OF RADIATION-DOMINATED FLOW 

Radiation pressure plays a crucial role in determining 
the dynamical structure of the accretion flows in bright 
X-ray pulsars. In order for the accreting gas to be bought 
to rest at the surface of the star by the radiation pres- 
sure gradient, the X-ray luminosity must satisfy the con- 
straint Lx ~ icrit, with the critical luminosity given by 
(Becker 1998; Basko & Sunyaev 1976) 



2.72 x lO 37 ^ 



ergs s 



(1) 



where Af* is the stellar mass, -R* is the stellar radius, 
ro is the radius of the cylindrical accretion column, er T 
is the Thomson cross section, and er and o~± denote 



2 
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the mean values of the electron scattering cross sections 
for photons propagating parallel or perpendicular to the 
magnetic field, respectively. The observation of many X- 
ray pulsars with Lx ~ 10 36 ~ 38 ergs s -1 implies the pres- 
ence of radiation-dominated shocks close to the stellar 
surfaces (White et al. 1983; White et al. 1995). 

The accretion flow is assumed to be cylindrical and 
highly supersonic in the upstream region. In this case the 
exact solution for the flow velocity v is given by (Becker 
1998; Basko & Sunyaev 1976) 



1 



— l+x/x st ' 



(2) 



The coordinate x increases in the direction of the flow, 
and v c is the flow velocity at the sonic point [x = 0), 
where the Mach number with respect to the radiation 
sound speed equals unity. The quantity x st represents the 
distance between the sonic point and the stellar surface, 
which is related to the other parameters via 



(3) 



In the far upstream region, x — > — oo and v — > (7/4) v c 
according to equation J5J • Since the incident (upstream) 
velocity is expected to be close to the free-fall velocity 
onto the stellar surface, we find that the sonic point ve- 
locity v c is given by 




4 (2GM, 



R,. 



1/2 



(4) 



3. CALCULATION OF THE RADIATION SPECTRUM 

In a steady-state, one-dimensional model, the Green's 
function / G representing the response to the injection 

of iVo photons per second, each with energy eo, from a 
source located at x = Xq satisfies the transport equation 

v 9f G = ^_i 21° + —( ° dJcA f a 
dx dx 3 de dx \ 3n e cr| dx 



to 



+ 



N 5(e - e )5(x - x ) 



(3v Q f G S(x - x ) 



(5) 



where e is the photon energy, n e is the electron num- 
ber density, vq = v(xq) is the flow velocity at the source 
location, and e 2 f G de gives the number density of pho- 
tons in the energy range between e and e + de. From 
left to right, the terms in equation (J3J) represent the 
comoving (advective) time derivative, first-order Fermi 
energization ( "bulk Comptonization" ) in the converging 
flow, spatial diffusion parallel to the column axis (i.e., 
in the cc-direction) , escape of radiation from the column, 
photon injection, and the absorption of radiation at the 
thermal mound, respectively. 

The quantity t csc appearing in equation JSJ) represents 
the mean escape timescale for photons to diffuse out 
through the walls of the accretion column. If the col- 
umn cross section is optically thick to electron scattering 
as expected, we find that t csc = ro/u>j_, where w± — c/t±_ 
is the diffusion velocity perpendicular to the n-axis, and 
t±(x) = n e (x)a±ro is the perpendicular optical thick- 
ness of the accretion column at location x. The dimen- 
sionless parameter /3 in equation JSJ controls the amount 



of absorption occurring at the thermal mound. In our 
application, f3 must be computed self-consistently since 
the mound is a blackbody surface that acts as both a 
source and a sink of radiation. We will use an energy 
balance argument to determine (5 below. 

Our primary goal is to calculate the radiation spec- 
trum / resulting from the injection of blackbody radia- 
tion at the thermal mound. Once the solution for the 
Green's function f G has been determined, we can easily 
compute / by convolving f G with a Planck distribution. 
If e 2 S(e)de denotes the number of photons injected into 
the accretion column per second in the energy range be- 
tween e and e + de from the blackbody surface at x = xq , 
we have (Rybicki & Lightman 1979) 



S(e) 



2ir 2 



1 



h 3 c 2 



(6) 



where Tq = T(xq) denotes the gas temperature at the 
mound "photosphere." The corresponding particular so- 
lution for / is then given by the integral (Becker 2003) 

f G {x ,x,e Q ,e 



f(e,x) 



e S(e )de 



(7) 



to N 

The accretion column in a bright X-ray pulsar is 
radiation-dominated, and therefore no net energy can be 
emitted or absorbed at the thermal mound. Instead, the 
absorption of radiation energy is almost perfectly bal- 
anced by the production of fresh blackbody emission. 
Most of the energy appearing in the emergent X-rays 
is transferred to the photons via collisions with the in- 
falling electrons while the radiation is diffusing through 
the column. Based on the energy balance requirement 
at the mound (x = xq), we can write the dimensionless 
absorption constant f3 in equation J5J as 

f °°e*S(e)de 



= 



(8) 



Trrlv U v U 

where Uq = U(xo) denotes the radiation energy density 
at the top of the thermal mound. 

Calculation of the radiation spectrum is simplified con- 
siderably if we work in terms of the new spatial variable 
y, defined by 

4 v(x) f 7 \- 1+x/x * 
y(x) = 1 - --U- = - , (9) 



which has the convenient property that y — > in the far 
upstream region (x — > — oo), and y — * 1 at the surface of 
the star [x — > x s t). Equation (J5J is separable in energy 
and space using the functions 



f\(e,v) 



g{\y) , 



(10) 



where A is the separation constant and the spatial func- 
tion g satisfies the differential equation 



2/(1 - y) 



1 



dy 2 



- 5y dg_ 
4 dy 

3(3v 5{y 



A 

4 + 
Vo) 



1 



4/y 



9 



7v n 



with solutions given by 



V) 



y <ya 
y>yo 



(ii) 



(12) 
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Fig. 1. — High-energy power-law sp ectr al index ao of the emitted 
photon number distribution <E> £ (eq. 1151 ') plotted as a function of 
the source location yo for the indicated values of the absorption 
parameter (3. Note the steepening of the radiation spectrum that 
occurs when fi is increased for a fixed value of yo, which reflects 
the decreasing residence time for the photons in the plasma due to 
the enhanced absorption. 



where y = y{x ) and B = Lp 1 (\ 7 y )/ip 2 {^,yo) to en- 
sure continuity of g at y = yo as required. The func- 
tions ipi and if2 are combinations of hypergeometric 
functions satisfying the appropriate upstream and down- 
stream boundary conditions (Becker & Wolff 2005). The 
system must also satisfy the eigenvalue relation 



dy 



<P2 



dipi 3f3v ip 1 (p 2 



dy 7v c y(l-y) 







(13) 



Equation i|13|l determines the eigenvalues A„ and the 
associated set of orthogonal eigenfunctions g n (y) = 
g(X n ,y). The first eigenvalue Ao is especially impor- 
tant because it characterizes the shape of the high-energy 
spectrum emerging from the column. 

Once the eigenvalues and eigenfunctions are deter- 
mined, the solution for the Green's function can be ex- 
pressed as the infinite series 



/ G (2/o,y,£o,e) 



E 

n=0 



en 



9n(y) , (14) 



where the expansion coefficients C n depend on the 
strength and location of the source. The phase-averaged 
spectrum emitted by the X-ray pulsar is computed by 
integrating the particular solution / over the height of 
the column to obtain 



„2 ,2 



t es lf(£,x) dx 



(15) 



where $ e (e) de gives the number of photons escaping from 
the column per unit time with energy between e and 
e + de, and / is computed using equation (JjJ. The asso- 
ciated high-energy photon spectral index, ao = Ao — 2, 
is plotted in Figure 1 as a function of the model param- 
eters (3 and yo- Note that a broad range of indices can 
be produced, although we must have ao > 2 in order to 
avoid an infinite photon energy density since the model 
considered here does not include a high-energy cutoff. 



4. ASTROPHYSICAL APPLICATIONS 

For given values of the stellar mass M* and the stellar 
radius our model has three free parameters, namely 
the column radius ro, the temperature at the top of the 
thermal mound To, and the accretion rate M. The first 
step in calculating the photon number spectrum $ c is to 
determine the value of xq (or, equivalently, yo) m terms 
of ro, To, and M, which fixes the location of the ther- 
mal mound. The velocity and density at the top of the 
mound, vo and po, respectively, are related via the con- 
tinuity equation M — irroVoPo- The density po can be 
calculated by setting the Rosseland mean of the free-free 
absorption optical thickness across the column equal to 
one, which yields for pure, fully-ionized hydrogen (Ry- 
bicki & Lightman 1979) 



Po = Jvo 1 = 4.05 x 10~ 12 T 7 



-12 rr.7/4 -1/2 



(16) 



in cgs units, where we have set the Gaunt factor equal 
to unity and introduced the conserved mass flux J = 
M/(7rrn). Once vo is computed, we can combine equa- 
tions (QJ and 10 to obtain 



2/o = l- v 



2GM, 



-1/2 



(17) 



which determines the location of the mound. For typ- 
ical X-ray pulsar parameters, we find that yo $ 1, and 
therefore the mound is situated very close to the surface 
of the star as expected. The final parameter that needs 
to be evaluated is 0. Using equations (@J and ©, along 
with equation (3.11) from Becker (1998), we find that the 
spatial variation of the radiation energy density is given 

by 



U(y) 



3My ( 2 GM, 



1/2 



Combining equations (JSJ and (|18fl . we find that 



= 



oT£ {2GM, 



3My v 



1 1 



-1/2 



(18) 



(19) 



Taken together, equations (QJ, (JTJJ, JT7I), and (JTHJl allow 
the determination of the four model parameters v c , vq, 
y , and (3 in terms of M„, r , T , and M. This closes 
the system and facilitates the calculation of the emergent 
photon number spectrum $ £ using equation 1151) . 

In Figure 2 we plot the theoretical predictions for the 
observed count-rate flux distribution, F e = $ e /(47r_D 2 ), 
along with the observed spectra for the X-ray pulsars 
4U 1258-61 (GX 304-1) and 4U 0352+30 (X Persei), re- 
spectively, where $ e is computed using equation H15JI 
and D is the distance to the source. The 4U 1258- 
61 data was reported by White et al. (1983), and the 
4U 0352+30 data is the result of XSPEC analysis of 
an archival RXTE observation taken in July 1998 and 
reported by Delgado-Marti et al. (2001). The observa- 
tional data represent the de-convolved (incident) spectra, 
and consequently there is some weak model dependence 
in these results. In both cases we set M* = 1.4 M© and 
i?» = 10 km. For 4U 1258-61, the model parameters used 
are r = 6 km, T = 7.3 x 10 6 K, M = 2.69 x lO^gs" 1 , 
and D = 2.5 kpc, which yields for the photon index 
ao = 2.04. For 4U 0352+30, we set r = 1.3 km, 
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1 2 5 10 20 



e (keV) 

Fig. 2. — Theoretical count-rate flux distribution F t (e) = 
<I> e /(47rD 2 ) computed using various amounts of interstellar absorp- 
tion are compared with the observed X-ray spectra of 4U 1258-61 
and 4U 0352+30. The column densities used are Nn = (solid 
lines); N K = 3xl0 21 cm~ 2 (dashed lines); and JV H = 9xl0 21 cm~ 2 
(dot- dashed lines). The other parameters for the theoretical models 
are given in the text. 



T = 9 x 10 6 K, M = 3.23 x 10 13 g s" 1 , and D = 0.35 kpc, 
which yields a.Q = 2.64. The distance estimates are taken 
from Negueruela (1998), and various amounts of inter- 
stellar absorption are included as indicated in the plots. 
The second source (X Per) is included due to its low lumi- 
nosity, which presents an interesting test for the model. 
Although the results presented here are not fits to the 
data, we note that the general shape of the spectra pre- 
dicted by the theory agrees fairly well with the obser- 
vations for both sources, including the turnover at low 
energies and the power law at higher energies. Several 
other sources yield similar agreement. In our model, the 
turnover at ~ 2keV is due to Planckian photons that 
escape from the accretion column without experiencing 



Basko, M. M., & Sunyaev, R. A. 1976, MNRAS, 175, 395 

Becker, P. A. 1998, ApJ, 498, 790 

Becker, P. A. 2003, MNRAS, 343, 215 

Becker, P. A., & Wolff, M. T. 2005, ApJ, submitted 

Burnard, D. J., Arons, J., & Klein, R. I. 1991, ApJ, 367, 575 

Coburn, W., et al. 2002, ApJ, 580, 394 

Davidson, K. 1973, Nature Phys. Sci., 246, 1 

Delgado-Marti', H., Levine, A. M., Pfahl, E., & Rappaport, S. A. 

2001, ApJ, 546, 455 
di Salvo, T., Burderi, L., Robba, N. R., & Guainazzi, M. 1998, 

ApJ, 509, 897 
Klein, J., et al. 1996, ApJ, 457, L85 
Meszaros, P., & Nagel, W. 1985a, ApJ, 298, 147 



many scatterings. This effect will tend to reduce the 
amount of absorption required to fit the observational 
data. 

5. CONCLUSIONS 

We have shown that a simplified model comprising a 
radiation-dominated accretion column with a blackbody 
source/sink at its base and a radiative shock is able to 
reproduce the power-law spectra observed in accretion- 
powered X-ray pulsars with a range of luminosities and 
power-law indices. The power-law plus blackbody form 
for pulsar spectra has previously been adopted in a 
purely ad hoc manner, but the new model described here 
provides a firm theoretical foundation for this general re- 
sult. Our work represents the first ab initio calculation 
of the X-ray spectrum associated with the physical ac- 
cretion scenario first suggested by Davidson (1973). 

Our goal in this Letter is to highlight the direct role 
of the accretion shock in producing the observed spec- 
tra, and therefore we have not included many other de- 
tails that are certainly relevant in some sources. For 
example, the observed spectrum of Her X-l has a very 
flat power-law shape up to a sharp exponential cutoff. 
The photon index in this case is a ~ 1, which is out- 
side the range allowed by our model since no high-energy 
cutoff is included here. In order to explain this shape, 
additional effects need to be included such as thermal 
Comptonization, cyclotron features, iron emission lines, 
and the possibility of an ultrasoft component produced 
by the accretion disk or the surrounding stellar surface. 
The inclusion of cyclotron cooling is particularly inter- 
esting since it connects the magnetic field with the flow 
dynamics by altering the pressure distribution, which can 
lead to further changes in the spectrum. We intend to 
explore these issues in future work. 
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